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r ~2' ABSTRACT 

We present a rigorous and practical way of constraining the Galactic potential based 
, on the phase-space information for many individual stars. Such an approach is needed 

to dynamically model the data from ongoing spectroscopic surveys of the Galaxy and 
in the future Gaia. This approach describes the orbit di stribution of stars by a family of 
parametrized distribution function (DF) proposed by McMillan & Binncy, which are 
based on actions. We find that these parametrized DFs are flexible enough to capture 
Jjj \ well the observed phase-space distributions of individual abundance-selected Galactic 

sub-populations of stars ('mono-abundance populations') for a disc-like gravitational 
potential, which enables independent dynamical constraints from each of the Galactic 
mono-abundance populations. 

We lay out a statistically rigorous way to constrain the Galactic potential pa- 
rameters by constructing the joint likelihood of potential and DF parameters, and 
. subsequently marginalizing over the DF parameters. This approach explicitly incor- 

porates the spatial selection function inherent to all Galactic surveys, and can account 
for the uncertainties of the individual position-velocity observations. 

On that basis, we study the precision of the parameters of the Galactic potential 
£N| ■ with that can be reached with various sample sizes and realistic spatial selection func- 

tions. By creating mock samples from the DF, we show that, even under a restrictive 
CN| ■ and realistic spatial selection function, one can recover the true potential parameters 

to a few per cent with sample sizes of a few thousands. The assumptions of axisym- 
metry, DFs that are smooth in the actions and no time-variation remain important 
limitations in the approach. 

x ■ 
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1 INTRODUCTION 

Constraining the Galactic potential optimally has been a 
very difficult yet important aspect of studies of the Milky 
Way and the Local Group. Understanding the various com- 
ponents of the Galactic potential, for instance separating 
the potential contributions from the baryonic disc and that 
from the dark matter halo is fundamental to understanding 
the history and formation of our own Galaxy. Furthermore, 
understanding the dark halo potential near the Sun is a cru- 
cial step to pin down the density and thus the scattering 
cross-section of the dark matter, which in turn is crucial in- 
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put for interpreting any dark -matter annihilation si gnal at 
the centre of Milky Way fe.g. lSu fc Finkbeineill20l3 l 

In the past few decades, constraining the Galactic po- 
tential h as mostly relied on the Je ans equation (for a sum- 
mary, see lBinnev fe T rcrnainc 2008) despite the known prob- 
lems of this approach. For instance, one only keeps the veloc- 
ity dispersion moments up to second order in the Jeans equa- 
tion. This predicts a Gaussian ellipsoid velocity distributions 
that does not capture the observed and expected asy mme- 
tries in the angular velocity (e.g. iFuchs et al.l l2009). 

One of the ways to alleviate these problems is through 
explicitly modelling the stellar population distribution func- 
tion (DF), either via integrals of motions or actions 
instead of just taking position-ve locity momen ts , as is 
the c ase of the Jeans equation jDehnenl 1 19991: iBinneyl 
l2010l : iBinnev fc McMillan! l201ll : iMcMillan fc Binnevl 12012 : 
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ISolwav et al.ll2012T ). For such modelling, it is crucial that rel- 
atively simple analytic distribution functions, e.g. based on 
the actions, are reasonably good approximations to the ac- 
tual orbit distributions for at least sub-populations of stars 
in the Galactic disc. How well such DF approximations work 
in practice is not yet clear. To explore this question, we 
consider in particular stellar sub-populations of a particu- 
lar age or abundance p attern (jFreeman fc Bland-Hawthornl 
l2002l ; lTing et al.ll2012l). 

Recently I Bow et al | (|2012al lbllch (hereafter, BO 12) 
showed, using G-dwarfs from the SEGUE survey (Sloan 
Extension for Ga lactic Understanding and Exploration; 
lYannv et al"1l2009l ), that mono-abundance, stars that have 
nearly the same [Fe/H] and [a/Fe], populations have very 
simple (a?, v) phase-space distribution properties: spatial 
density distributions that vary exponentially, both in the 
radial and vertical direction; and furthermore, each mono- 
abundance population shows isothermal velocity dispersion 
in the v ertical direction an d exponentially in the radial di- 
rection (|Bovv et al.l l2012ch . This simplicity in [x, «)-space 
gives hope that one could describe mono-abundance pop- 
ulations with simple action-based DFs which is one of the 
purposes of this study. 

Since stars are the most abundant tracers of the Galac- 
tic potential, simple action-based DFs will naturally provide 
a statistical and rigorous method to constrain the Galactic 
potential because the conversion of position-velocity vari- 
able (configuration space) to the actions-angl es variable de- 



pend s on parameters of the po t ential (see iBinney 201(j 
20 111 : iBinnev fc McMillan! l201lt iMcMillanl (201 it (BinnevI 
20121 ; [McMilla nfe BinnevI |2012| . for generic ground-work 
on the topic). The main idea is to determine the likeli- 
hood of the observational data, given a joint set of param- 
eters for both the DF and the gravitational potential; sub- 
sequent marginalization over the DF parameters provides 
then a rigorously derived constraint on the potential. Such 
an approach appears to be a precondition to fully exploit 
the dynamical information content of large Galactic sur- 
veys that will provide us spatial and velocity distribution, 
along with elemental abundances, including APOGEE (The 
Apache Point Observ atory Galactic Evolution Experiment; 
lEisenstein et all 1201 ll ) , GAL AH (Ga l actic A rchaeology sur- 
yey with HERMES ; iFreemanl |2010| . |2012j ) and ESO-Gaia 
l|Gilmore et alj|2012h . 

In the context of dynamical modelling the Galactic disc 
with action-based DFs, the goals of the paper are three- 
fold: explore how well 'mono-abundance populations' can be 
approximated by simple action-based analytic DFs; lay out 
a formalism that provides constraints from a set of discrete 
stellar positions and velocities on the gravitational potential, 
after marginalizing over the DF; forecast what constraints 
on the shape of the potential we can expect given existing 
sample sizes. 

This paper is organized as follow: in Section [5] we dis- 
cuss the choice of parametrization of the DF and summarize 
the adiabatic approximation that we assume to calculate the 
actions variable efficiently. In Section [3J we will show that 
our choice of DF exhibits all the basic phase-space distri- 
bution properties of mono-abundance stars, with the only 
caveat that disc potential has to be included in order for the 
vertical spatial profile to fit. In Section |4j we show that by 
studying the likelihood of the DF-potential parameters, one 



can recover the potential parameters, even under a restric- 
tive and realistic spatial selection function. 

Throughout this study, we denote cylindrical coordinate 
to be x = (R, 4>, z). We assume velocities are measured in the 
inertial Galactocentric frame with axisymmetric potential, 
$ = $(R,z). 



2 CHOICE OF DISTRIBUTION FUNCTION 



We follow the DF advo cated by IBinnevI (|201Ch and 
IMcMillan fc BinnevI (|2012l ) closely, but note that we rear- 
range some of the terms in the DF to facilitate physical in- 
terpretation of the DF. For circular orbits in the equatorial 
Galactic plane, the DF has to be of the form: 



f(Jn, L z , J z ) = f(Lz)8(Jn — 0)S(J z — 0) 



(1) 



where f(L z ) determines the angular momentum distribution 
for a infinitely cold planar disc, Jr and J z are the radial and 
vertical action, respectively. 

In the cold disc limit, the circular radius R c coincides 
with the observed radius, where 



L Z = R C - V C (R C ) 

and V c is the circular velocity, 



R 



9$ 
dR 



(2) 



(3) 



Thus we have for an exponential disc, up to a multi 
plicative normalization constant, 

R, 



f(L z )dL z 



exp 



( _Rc 
V h R 



RcdR c 



(4) 



The DF parameter fin determines the radial spatial distri- 
bution. In a more general setting, we relax the Jr = J z — 
constraint and propose the DF to be: 



f{Jn,L»J,) = f(L z ) 



where 

Cr(L z ) = oh exp 



D R (L Z 



a z exp 



k{L z ) 
Cr(L z , 

D R (L Z 

Rg{L z ) 

Rc(L z ) 



exp 



■ exp 



k(L z )Jr 
Cr(L z ) 

v{L z ).] z 
D R (L Z ) 



(5) 



(6) 



(7) 



where R c is the circular radius given angular momentum 
L z , assuming equatorial in-plane movement, k(L z ) and 
v(L z ) are the radial and vertical frequency, respectively, 
under epicycle approximation. The DF parameter h a deter- 
mines the radial exponential decay of the velocity dispersion, 
whereas <tr and a z control the total radial and vertical dis- 
persion, respectively. 

The n and v terms are necessary (e.g. lBinnevll2010h for 
two reasons: 

(i) A finite increment in energy E can lead to an infinite 
increment of actions in the unbound case. One would love 
to have a term that will tend to zero at large L z and couple 
this term with Jr. One of the possible choices is the radial 
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frequency k(R c (L z )). At large R c , the effective potential is 
less concave at the minimum point and therefore n(R c (L z )) 
tends to zero. Similarly for the vertical frequency. 

(ii) Qualitatively, in R c region where k changes more 
drastically, the increment of Jr is also more drastic, and 
therefore the scale length of Jr has to decrease in propor- 
tion to compensate this effect. Similarly for the vertical os- 
cillation. 



2.1 Adiabatic approximation 

Importantly, we only measure (-R, z, Vr, v$, v z ) but not the 
actions. In order to study actions-based DF, we employ the 
adiabatic approximation in order to calculate the radial and 

vert ical actions efficiently . 

iBinnev fc McMillan! (|201l| j and ISolwav et all (|2012r i 
showed that the vertical action J z calculated in a fixed ra- 
dius R approximation is almost conserved, even with radial 
migration, and interferences of spiral(s) or bar. More pre- 
cisely, 



■Jz 



dz z(z\R) 



(8) 



This reduces the calculation to one-dimensional integra- 
tion that can be numerically calculated easily. Furthermore, 
since $(z\R) = <&(—z\R), to calculate the approximated ver- 
tical action, we have: 



Jz 



dz J%~^ +2^(z lalt .\R 1Dlt .) 



(9) 



where the subscript 'init.' stands for the observed value in 
practice. 

Similarly, for the radial action, we can consider the ra- 
dial component at the equatorial plane, ignoring the vertical 
motion. Under this assumption, 



Jr = * dR R(R\z = 0) 



We deduce 



Jr = - 



dR 



RL, + 3 + 2 (* (iitalt ' 0) " m ' 0) ) 

where 7? apo . and iZ pori . are the apo- and peri-centric radius, 
respectively. 



2.2 Variation of actions under adiabatic 
approximat ion 

Fig [T] and Table [T] illustrate the variation of the approxi- 
mated actions in the course of a trajectory. We assume a 
Miyamoto-Nagai potential, 

GM 



(12) 



R 2 + (a+ Vz 2 + b 2 ) 2 



with parameter GM — 7.5 x 10 5 km 2 /s 2 kpc, a — 5 kpc, 
6=1 kpc. We choose these parameters because we have 
almost a flat rotation curve at Galactocentric radius 5 kpc < 
R < 12 kpc and circular velocity V c (Rq) — 220 km/s at the 
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Figure 1. Approximate action calculation in a Miyamoto-Nagai 
potential: the top panel shows the circular velocity of a Miyamoto- 
Nagai potential with GM = 7.5 X 10 s km 2 /s 2 kpc, a = 5 kpc, 
6 = 1 kpc, the vertical dashed line marks approximately the so- 
lar Galactocentric radius, the horizontal dashed line indicates the 
common accepted circular velocity V c = 220 km/s at solar radius 
and the shaded region shows the region where the circular velocity 
is approximately constant in this choice of potential; the middle 
and bottom panel show the variation of actions calculated using 
adiabatic approximation in the course of trajectory. We plot the 
case of -Ri n it. = 6 in Table [TJ The results show that actions cal- 
culated with adiabatic approximation arc almost conserved along 
the orbit with less than 5 per cent variation. 
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Table 1. Variation of actions variable under adiabatic approx- 
imation with slight perturbation to the circular orbit under 
Miyamoto-Nagai potential that gives a flat rotation curve. 



Riait. 






#init. 




Aj 2 


kpc 


km/s kpc 


km/s 


km/s 


per cent 


per cent 


3. 


450. 


90. 


60. 


4 


3 


6. 


1250. 


60. 


10. 


5 


3 


9. 


1970. 


30. 


20. 


4 


2 



solar radius Rq = 8 kpc as shown in the top panel of FigQ] 

We apply slight perturbation on circular orbits that 
mimic the radial velocity dispersion of 40 km/s at solar 
Galactocent ric radius, radial t o vertical velocity dispersion 
ratio of y/2 jFuchs et all 120091 ) and velocity scale length of 
3.5 kpc (B012). We define the variation of action Aj to 
be the one-sigma over the mean value of the action for each 
time-step. The results of A j are listed in Table[T] The results 
show that the approximated actions are almost conserved in 
all cases. 

We also examine potential due to double exponential 
disc (i.e. with density decaying both radially and vertically), 
and also flattened isothermal potential: 



* = V c 2 (R@)ln 



S1 + ? 



(13) 



The results are qualitatively similar, with the vari ation of 
actions less than 5 per cent in all cases. lBinnevl(|2012h studies 
a more precise method of calculating the actions, but as we 
will show in the results below, it is not necessary for our 
study purpose. 



3 MATCHING THE PHASE-SPACE 

DISTRIBUTION OF MONO-ABUNDANCES 
POPULATIONS 

In this section, we will show that the DF in equation (J5J) is 
a good representation of the observed position-velocity dis- 
tribution of mono-abundance sub-components of the Milky 
Way disc, namely that they exhibits an exponential spatial 
density, both radially and vertically, and an isothermal ve- 
locity dispersion in the vertical direction and exponentially 
decaying in the radial direction. 



3.1 Dealing with Jacobian change of measure 

In order to have the conversion of DF between canonical 
position-velocity variable and the actions variable, one has 
to take care of the Jacobian term carefully. Note that 



d 3 xd 3 pf(x,p) = / d 3 Jd 3 8f(J) 



(14) 



where p is the canonical momentum associated to x. Since 
(x,p), (J, 9) axe canonical variables and measure is invari- 
ance under canonical transformation, we have 



Furthermore, in cylindrical coordinate, (pR,p<t,,p z 
(vr, Rv$, v z ), we deduce 



f{x,v) = Rf(j(x,p{x,v))^j 



(16) 



In summary, for an axisymmetric potential, in order 
to convert actions-based DF to (x, «)-space distribution, it 
suffices to multiply the Jacobian (= R) term. 



3.2 Fitting the phase-space distribution via 
maximum-likelihood 

We create mock data stellar tracers with (R, z, vr, v$, v z ) 
attributes, such that it has density scale length about 2.5 
kpc, dispersion scale length about 3.5 kpc with Gaussian 
velocity dispersion of or — 40 km/s at solar Galactocen- 
tric radius, <tr/o z = \[2, orIo^ = 1.5 l)Fuchs et al.ll2009l . 
B012). We assume that the mean velocity of the Gaussian 
for each spatial point 



vW(R, z) = v;{R, z) = 



Vj(R, z) = vj(R) = V C {R) - A 



2V C (R) 

where (see I Bow fe Tremainell2"oi2l . for details), 



A 



(— + — 

\ h R h a 



x R-0.5 



(17) 



(18) 



(19) 



We assume that all of the mixed moments (o-fj^, tr Rz 
and o-0 Z ) vanish. Thus, we assume that both the vertex de- 
viation and the tilt of the velocity ellipsoid are zero, even 
though these were not measured by B012. The tilt of the 
velocity ellipsoid in particular is not e xpected to be zero a t 
heights > 1 kpc above th e plane (e.g. [Siebert et al.l l2008h . 
However, as discussed by IBinnev fc McMillan! ( 201lf ). the 
adiabatic approximation is unable to capture a non-zero tilt 
as it assumes that the radial and vertical motions are in- 
dependent. Thus, even though the quasi-isothermal DF of 
equation ([5]) has a non-zero tilt when used with correctly 
calculated actions (for exa mple, using the torus machinery; 
IBinnev fc McMillan! l20l"lh , it does not when using the adi- 
abatic approximation. For this reason, we assume that the 
tilt of the velocity ellipsoid is zero in what follows. 

We search for the best-fitting DF parameters by maxi- 
mum likelihood. For this section, we fix a particular choice 
of Miyamoto-Nagai potential parameters: GM — 7.5 x 
10 s km 2 /s 2 kpc, a — 5 kpc, b — 1 kpc. We examine var- 
ious choice of the potential parameters. It does not change 
the results qualitatively. We emphasize that the mock trac- 
ers population mimics a mono-abundance population. Note 
that, the probability that we observe (x a ,v a ), according to 
our DF model with DF-parameters A DF IS 



Pr(f Q , Wa|$, Adf) = f(x a ,V a \$, Adf) 



(20) 



Given this probability, we can define the log-likelihood 
to be the sum over N a mock data, and search for the best 
DF-parameters by maximizing the log- likelihood. 



In £df(Adf) = ^ ln[Pr(f a ,v a \$, A DF 



(21) 



f(x,p) = f(j(x,p)) 



(15) 



We obtain the best-fitting parameters by performing 
a nested-grid search on the multi-dimensional Adf space 
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Figure 2. How well can the DF family from equation {5} approximate the observed properties of Galactic mono-abundance stellar 
sub-populations (B012) in a Miyamoto-Nagai potential? Those mono-abundance populations have an exponential spatial density, both 
radially and vertically, and an isothermal velocity dispersion in the vertical direction and exponentially decaying in the radial direction. 
To create mock data, we assume a radial velocity dispersion about 40 km/s at solar Galactocentric radius and scale height h z = 0.18 
kpc for the mock data. The top row shows the radial (left) and vertical (right) density distributions in (R,z) £ [0,12] X [—0.5,0.5] kpc 2 
and (R, z) S [7, 9] X [0, 1.3] kpc 2 , respectively. The dashed red line shows the reconstructed profile from the best-fitting DF and the grey 
solid histogram shows the mock data. The second row shows the distribution of vr as a function of R, with z 6 [—0.5,0.5] kpc (left 
panel) and as a function of z with R S [7,9] kpc (middle panel), and also the histogram in (R, z) S [7,9] X [—0.5,0.5] kpc 2 . The grey 
shaded region shows the reconstructed profile from the best-fitting DF with the white solid line to be median value, whereas the black 
solid filled circle and their vertical error bars show the results of the mock data. The third and fourth row show the distribution of v z 
and i>0, respectively. The DF family fits the tracer population very well for a disc potential. 
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Figure 3. Similar to the top right panel of Fig. [2] but for isother- 
mal potential with the axis ratio q = 1 and circular velocity 
V c = 220 km/s. We also tried various choices of q, but the fits are 
equally unsatisfactory: the DF family in equation J5jl cannot pro- 
duce a vertically exponential tracer density profile, if stars orbit 
in a flattened isothermal sphere. 



and finding the set of parameters that maximizes the like- 
lihood. After obtaining the best-fitting parameters, we cal- 
culate from f(j(x,v)) the stellar profile in configuration 
space (x, v) and compare it to the one of the mock data 
to check whether our assumption on the functional form of 
the DF represents well the mock data. We consider a grid 
of position-velocity configuration space (R, z, vr, v$,v z ) and 
calculate the DF in configuration space using equation (116[) . 
We then calculate the stellar profile of each position- velocity 
component by marginalizing over the other spatial-velocity 
components. 

Note that, given the functional form of the DF in equa- 
tion (JSJ, the predicted vertical profile is fixed by the verti- 
cal velocity dispersion and the potential, but independent of 
the scale height of the mock data. We find that the vertical 
profile is predicted satisfactorily with the Miyamoto-Nagai 
potential (as shown in Fig. [2| but not with the isothermal 
potential (as shown in Fig. 

This can be explained by the following. In the case that 
R^> z, the vertical restoring force of the Miyamoto-Nagai 
potential and thus z is about constant. This implies that the 
vertical action is proportional to z, and the vertical profile in 
ansatz equation ([5]) goes down exponentially. On the other 
hand, the vertical restoring force of the isothermal potential 
and thus i is proportional to z. This implies that the vertical 
action are proportional to z 2 , and the vertical profile will 
goes down too drastically [proportional to exp(— z 2 )]. 

We conclude that the quasi-isothermal DF of equa- 
tion ([5| provides a good representation of the DF of mono- 
abundance sub-populations of the Milky Way disc. More 
importantly, the action-based DFs predict a low-v^ tail. As 
discussed in Section [[] this co incides better with the obser- 
vation (e.g. iFuchs et al.1 12009). In the following section, we 
will generate mock samples by rejection sampling directly 
from the DF. 



4 CONSTRAINING THE GALACTIC 
POTENTIAL 

In this section, we will show that, given a functional form 
of the potential and the DF, we can recover the optimal 
potential parameters by marginalizing over the DF parame- 
ters. Note that in practice, we only observe a very small part 
of the Galaxy. For simplicity, we assume that the observed 
volume is a cylinder around the Galactic center with a finite 
width and height. Mathematically, we have to include the 
selection function 



1, if \R - 8| < 1 kpc and \z\ < 2.5 kpc 
0, otherwise 



(22) 



With this selection function in place, given any param- 
eters of the DF and the potential, A = (A DF , A pot ential ) , We 
carefully normalize the DF through Monte-Carlo integra- 
tion, more precisely 



where 



f(JW 



— f d s x d : 
2tt J 



vR-P aelect (R,z)f(J(x,v)\X) 



(23) 



(24) 



To our best knowledge, this is the first attempt that 
has never before been implemented to study the joint dis- 
tribution of DF and potential parameters in the context of 
constraining the Galactic potential. 

First, we consider the two-dimensional in-plane 
Miyamoto-Nagai potential which contains two potential pa- 
rameters: GM and a. We consider the case where ha — 
2.5 kpc, ha = 4 kpc and or = 12000 (km/s) 3 in equa- 
tion ©, which roughly mimics the observed velocity distri- 
butions from our calculation in Section [3] 

We define the likelihood of the potential parameters, 
£ = £(A P otc„tiai) by marginalizing over A DF - We estimate the 
one-sigma significance of the best-fitting potential parame- 
ter by finding the demarcation line at 2 ln(£/£ max ) = 2.7 
and two-sigma significance at 21n(£/£ max ) = 4.6. We create 
data with GM = 7.5 x 10 5 km 2 /s 2 kpc, a = 5 kpc, 6=1 
kpc in the observed volume R = 7 — 9 kpc and \z\ < 2.5 kpc. 

We reiterate the process for various sample sizes for 
the mock data, iV traoBra for the same observed volume. 
As shown in Fig [3] the one-sigma uncertainty goes down 
as l/V iVtracers- We use the results with the sample size 
iVtracers = 30 000 to be our reference, because the uncertainty 
of the potential parameters is ^ 0.1 per cent in this case. We 
show that this reference value coincides with the true value, 
indicating that, assuming a potential model, we can recover 
the true potential parameter with observation restricted to 
a small volume. We show that we can find the best-fitting 
potential parameters within one- to two-sigma significance 
from the true value regardless of the sample size. 

We extend this results to the three-dimensional case 
by considering the flattened isothermal potential with two 
parameters [equation (|13[) ] : V c and q. We consider the warm- 
disc case where ha = 2.5 kpc, h a — 3.5 kpc, or — 
12000 (km/s) 3 , o z = 8000 (km/s) 3 for DF in equation ©. 
We create the mock data from the potential with V c = 
220 km/s and q = 0.8. As shown in Fig[5j one can recover 
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Figure 4. Likelihoods of the data given the parameters of a Miyamoto-Nagai potential, assuming that the iVtracara observational 
constraints are available within R = 7 — 9 kpc and \z\ < 2.5 kpc. We only plot the likelihood region that encloses 95 per cent of the 
total probability. The left panel shows the results of mock data points iV trac< , rB = 1 000, the middle panel Af tracer3 = 3 000 and the right 
panel Af tracora = 30 000: the likelihood contours tighten as l/y/N tr ^ cers around the correct input parameters for the mock data. The data 
constrain a higly degenerate combination of GM and a. 
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Figure 5. Similar to Fig [4] but for a flattened isothermal potential. Here the estimates for the potential flattening and the circular 
velocity are hardly correlated. 



both the circular velocity V c and the axis ratio q successfully 
with our proposed method. 

Note that, in most cases, without a Galactic dust-map, 
observations in the equatorial Galactic disc are less reliable. 
Therefore, we also tested the more realistic scenario, where 
the observed volume is restricted to the annulus height of 
1 < z < 2.5 kpc without the mid-plane. The results show 
that we can about equally well constrain the axis ratio q, 
given the same sample size. This could be probably due to 
the fact that our three-dimensional isothermal potential is 
too restrictive, therefore the axis ratio is not too sensitive 
to the vertical observed volume restriction. However, the V c 
uncertainty is larger, indicating the adiabatic approximation 
is not as good away the from Galactic plane. 



5 CONCLUSIONS 

We demonstrate a way to constrain the Galactic potential 
from mono-abundance stellar populations. We show that 
the phase-space distribution properties of mono-abundances 
stars can be fitted very well with simple action-based distri- 



bution function. Assuming that the vertical and radial ve- 
locity dispersion scale lengths are the same, these properties 
are fully-determined by only four DF parameters. 

We further show that, assuming that the proposed DF 
is the true representation of the mono- abundance popula- 
tions and certain parametrization of the Galactic potential, 
we devise a statistical and rigorous method to measure the 
Galactic potential parameters. This is achieved by calculat- 
ing the likelihood of observational data, given a joint set of 
parameters for both the DF and the gravitational potential 
and subsequently marginalizing over the DF parameters. 

We create mock data sample from the distribution func- 
tion with various potentials, including the Miyamoto-Nagai 
in-plane potential and a flattened isothermal potential. We 
show that even with a mono-abundance sample of the size of 
a few thousands, we can pin down the potential parameters 
to a few per cent. 

In this study, we assume that the measurement is per- 
fect without uncertainty in the position and velocity mea- 
surements. The inclusion of measurement uncertainty should 
be conceptually straightforward. We will explore this aspect 
in a future paper. 
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Figure 6. Similar to the middle panel of Fig. [5] but with the ver- 
tical range of available tracer constraints restricted to 1 < z < 2.5 
kpc- the mid-plane cut out: this is a more realistic selection func- 
tion since in most cases, without a Galactic dust-map, observa- 
tions in the equatorial Galactic disc are less reliable. 
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Although a more realistic potential, for instance a com- 
bination of halo, bulge and disc potential will contain more 
than two parameters and hence post a challenge to the esti- 
mation, different mono-abundance stars are tracing the same 
Galactic potential and therefore their joint constraints with 
the method we propose is still very promising. 
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